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Abstract 



The Davydov model of energy transfer in molecular chains is reconsidered 
assuming the distance dependence of the exciton hopping term. New equa- 
tions of motion for phonons and excitons are derived within the coherent state 
approximation. Solving these nonlinear equations result in the existence of 
Davydov-like solitons. In the case of a dilatational soliton, the amplitude and 
width is decreased as a results of the mechanism introduced here and above 
a critical coupling strength our equations do not allow for localized solutions. 

For compressional solitons, stability is increased. 
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I. INTRODUCTION 



Davydov |I[ formulated a model intended to explain an almost lossless energy transfer 
in quasi-one-dimensional biomolecular chains such as peptides, proteins and the DNA. The 
key element in the Davydov model is the coupling between excitons and phonons of the 
biomolecular chain. The starting point is typically the Hamiltonian given by 



2 ^ 

^ -n 



+ £ [DAlA n - J(Al +1 A n + A n+1 A{)] (1.1) 
j n 

where the position of the n-th unit is 

x n = dn + u n (1.2) 

with d denoting the equilibrium separation and u n a displacement from equilibrium. Here, 
p n is the associated momentum, W the elastic modulus and M the mass of the oscillating 
unit. The operators A\ and A n create and annihilate an exciton at site n, respectively. 

The constant D describes the energy of a single exciton on a given site while J refers 
to hopping of excitons between neighboring sites. Since both constants originate from the 
dipolar nature of the chain's units it is obvious that they must be distance dependent. Hence, 
denoting the relative distances between two neighboring groups as 

r n,n+l — \x n +l ~ %n\ — d + Mn+1 — U n (1-3) 

and 

^n,n— 1 \%n %n— 1 1 d -\- U n 1l n _\ (1-4) 

we assume that 

D = Q + L>(r n , n+1 ) + D(r n , n _i) (1.5) 

and analogously that 

J = J{r n ,n+i) + J{r n ,n-i)- (1-6) 



In most of the literature on the subject only the dependence of D on the distance has been 
explicitly introduced. It is quite obvious, however, that J should be inversely proportional to 
the third power of the distance (dipole-dipole interaction) and this for small displacements 
gives an approximately linear dependence of J on the magnitude of the displacement. Hence, 
we can postulate that: 

D(r rhn+1 ) = D + xiiun+i - u n ) (1.7) 
£>(r njn _i) = D + xi(u n ~ t*n-i) (1-8) 

J{r n ,n+l) = Jo + X2<>n+1 - U n ) (1.9) 
J{r n ,n-l) = Jo + X2(u n - Un-x) (1.10) 

For a treatment of similar terms in a somewhat different spirit we refer the reader to Ref. 

First of all, let us note that xi is responsible for the existence of a localized solution 
called a Davydov soliton, and this fact is irrespective of the sign of Xi 0- F° r Xi > 0; the 
molecular chain is locally compressed (compressional soliton) and for xi < 0, the soliton is of 
dilatational type (locally dilatated). This is not an uncommon situation in soliton bearing 
nonlinear lattices |§. It is still an open question as to the value and the sign of Xi ( see 
for instance Ref. [Q]). The constant Xi is clearly negative since the interaction diminishes 
with distance. Scott || had estimated for peptide chains xi from the approximate relation 
^ 2 = 3J/d using J = 7.8 cm -1 and d = 4.5 A as X2 = lpN. He also gave an estimate 
of xi — 34piV. These values will be used in subsequent numerical calculations. We will 
show in the sections that follow that keeping both xi an d Xi m the Hamiltonian results in 
equations of motion which differ from the standard Davydov system. It will be of interest to 
examine their solutions and find out whether solitons can be formed and what affects their 
properties. 



3 



II. THE MODEL HAMILTONIAN 



Based on the results of the previous section we rewrite the Davydov Hamiltonian's exciton 
part as 

H exc = + 2A)) + - «n-l)] (2-1) 

n 

~J2 A l i( J + X2(u n - M n _!))A n _! + (J + X2(u n+1 - U n ))A n+l ] 

n 

where the first term is responsible for on-site excitations and the second for hopping. 

In the next stage we find the continuum limit of both the excitonic part H exc given above 
and the phonon contribution H p h. To make the continuum limit we take: 



a, s * a, x , dA , Id 2 A l2 

A n ^A(x); A n±1 = A(x)± — d+-—d 2 + ... 

/ \ ~ t \ , 9u , 1 (9 2 M , 2 

m„ -> M n± i = U{X) ± fa d+ 2d^2 + " ' 



(2.2) 
(2.3) 



and define 



t/ \ t du . 
J[x) = J + a 

(9m 

L>(x) = Do + xi-^ d. 



(2.4) 
(2.5) 



As a result, the phonon part of the Hamiltonian becomes: 

2/ n /o \ 2 



- I dx 

2 



p (#) 
M 



The exciton contribution to the Davydov Hamiltonian is found as 



Hp 



dx 



QA^x)A(x) + 2 Xl d(^\ A\x)A(x) 



-JdxU + X2 d[f x 



2A\x)A(x) 



dtf dA 

dx dx 



where we have denoted 



(2.6) 



(2.7) 



Q = fl + 2D 



(2.8) 



and the hopping part has been obtained from a symmetric form given by 



H hop = -E J (AlA n+1 + A[ +1 A n ) (2.9) 

n 

which is explicitly hermitian. 

Our next task is to derive equations of motion for the coupled phonon field u(x) and 
the exciton field A(x). This can be done using several approaches, e.g. through functional 
minimization, coherent state Ansatz, etc. In the following section we derive these equations 
by treating phonons classically and taking a special Ansatz for excitons. 



III. DERIVING THE EQUATIONS OF MOTION 

In order to derive equations of motion we assume the phonons to be virtually classical so 
that Hamilton's equations can be applied to u(x) and p(x). The exciton part is approximated 
by the expectation value of H exc in the trial ground state 



dx a(x)A\x)\0) 



(3.1) 



where a(x) is a complex probability amplitude of an exciton wave and |0)is the vacuum. 
Thus, the effective Hamiltonian used is 



H eff = H ph + (^j\H e 



(3.2) 



or explicitly: 



Heff = ~ I dx 



+ dx 



M \dx 



'-2J + n) + 2(xi-X2)d— 



(3.3) 



\a(x)\ 2 + (J + X2 d^] 



da 



dx 



d 2 



We now derive the equations of motion for the phonon field as 

--{H eff ,p}- — 



and 



2(xi-X2)d\a\ 2 + X 2d 3 



da 


2" 


dx 





(3.4) 



(3.5) 



where the curly brackets above indicate Poisson commutators. Differentiating eq. ( |3.4|) with 
respect to time and utilizing eq. (|3.5|) we obtain 



d 2 



d 2 



dt 2 V °dx 2 



u(x, t) 



l_d_ 

Mdx 



2(xi -X2)d\a\ 2 + X 2 d 3 



da 


2" 


dx 





(3.6) 



which is a wave equation for the displacement field with an additional term arising due to 



the coupling with excitons. The sound velocity v is defined as v = d^JW/M. 

We treat excitons quantum mechanically and hence use the Schrodinger equation 

d , , 



ih- 



'dV 



as an equation of motion. Here, H exc is taken in the form given in eq. 
verify that, as a consequence, we obtain 



ih 



da 
dt 



(n-2J +) + 2( Xl - X 2)d 



' du 
dx 



a\x) 



d_ 

dx 



(J + X2 d 



du /da" 
dx \ dx 



(3.7) 
It is easy to 

d 2 . (3.8) 



We note that the two coupled equations, eqs. ( pT6| ) and ( |3.8j ) resemble the Davydov set of 
equations. In fact, as \2 the result agrees identically with Davydov. 

Our next task is to solve these two equations with a specific objective of finding soliton 
solutions. 



IV. SOLVING THE EQUATIONS OF MOTION 

In the first step we postulate that 

u(x, t) — u(x — vt) (4.1) 

i.e. that the (real) displacement field is a travelling subsonic wave (v < vq) and 

a(x, t) = exp[— iEt + ia(x — vt)]a(x — vt) (4-2) 

which means that the complex exciton field travels with the same velocity v but also has 
internal oscillations. The energy of the soliton is not E, since in our Ansatz, we have other 
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time dependent pieces. A more complete form is easy to obtain 0. As a result, the two 
equations become ordinary differential equations in the new coordinate £ = x — vt. Thus, 



{v 2 - vl)u = — [2(xi - xa) da 2 + X 2 d 3 {a') 2 



(4.3) 



and 



-E + Q- 2 J 



Ad 2 J ( 



+ 2d(xi - X2)u 



a-d 2 



(Jo + X2 du)a 



(4.4) 



where the prime denotes differentiation with respect to £. Note that as X2 — ► 0, eqs. ( |4.3| ) 
and ( |4.4| ) reduce to the original system of Davydov equations . 

To proceed we note that the second term in the square bracket of eq. ( f4.3| ) is much smaller 
than the first one. This is on the account of the fact that for large solitons, a 2 3> (a) 2 - We 
can then integrate eq. ( [4.3| ) once to yield 



' ^ 2(xi -X2)d 2 
U ~ M(v 2 - vl) a ' 



(4.5) 



This can then be substituted into eq. (4.4) to give 

Aa 



d a — 



+ —(a ) 2 = 0. 



J + 7a 2 J + 7a 2 J + 7« 2 
where, to simplify notation, we introduced the following symbols 

2X2 (Xi ~X2)d 2 



(4.6) 
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M[v 2 - vl) 



4(xi-X 2 ) 2 ^ 2 

M{v 2 -v 2 ) ^ Uld 2 



■1 



(4.7) 
(4.8) 



and 



A = -E + Q - 2 J 



4J rf 2 



> 0. 



(4.9) 



Eq. (4.6) has two free parameters, E and v (the velocity of the soliton). The parameter 
E is fixed by the requirement that the excitonic wave function be normalized, ie. 



— / a 2 (z) dz = 1. 
d 



(4.10) 



Using the transformation 



we obtain from eq. (4.6) 



where 



and 



y(a) = % (4-11) 



yy' + f(a)y 2 + h(a) = (4.12) 



/(«) - (4.13) 



, . . Aa + Ta 3 , 

J + 7« 



Eq. ( |4 . 1 2|) is a Bernoulli equation which can be completely integrated and its general solution 
is 

y 2 = exp(-J)(ai - 2 F h(a) exp(J) da) (4.15) 

7 an 



where 



In our case 



1 = 2 F f(a)da. (4.16) 

J do 



1 = 21n(J + 7 a 2 ) (4.17) 



and eq. (|4.15|) yields 



'da\ ai + J Aa 2 + i(A 7 + J r)a 4 + i 7 ra 6 



= H(a). (4.18) 



Vdey (^ + 7a 2 ) 2 

This final equation can be integrated numerically in terms of periodic and localized solutions. 
Note that for a% 7^ only singular or periodic solutions can be found. Before we present a 
complete analysis of the various solutions we note that in the limit of X2 —* we find the 
standard elliptic equation f7j, i.e. 

8 



since 7^0. Here, elliptic waves can be found in a standard way, and for c = we find the 
Davydov soliton as 



a(fl = J-^sech(J^). (4.20) 



r \v j, 







V. ANALYSIS OF SOLUTIONS 



In solving eq. (|4.18|), we assumed the values of model parameters consistent with those 



used for peptide chains These values are: 

W = 40 N/m 
M = 5.7- 10~ 25 kg 
d = 4.5-10- 10 m (5.1) 

n ~ o.2- 10" 19 j 

J = 1-55 ■ 10~ 22 J 

For illustration purposes, we have chosen to investigate slow solitons (v/vq ~ 0.02). Physical 
values for \2 are on the order of —1 pN, and we will look at values within this range. Note 
that X2 must be negative since dipole-dipole interaction decreases with separation. The 
value of Xi is still an open question. Experiments seem to indicate that Xi — 30 — 60 pN. 
This would correspond to a compressional soliton. On the other hand, theoretical estimates 
range from —60 to + 25 pN raising the possibility of dilatational solitons (for xi < 0). 
Therefore, we will investigate both positive and negative values for xi- 

In the case of compressional solitons (xi > 0), qualitatively nothing special happens. 
Including the effect of X2 simply enhances the stability of the soliton. This can be understood 
as follows. With xi < 0, we see from eq. (4.7) that 7 > which immediately implies that 
the denominator of H(a) given in eq. (4.18) is increased. This means that in the compressed 
chain, the effective Jo is larger. 
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If xi < 0, it follows that 7 < 0. In this case, the soliton is dilated and the dipole-dipole 
interaction is weakened which can lead to a destabilization of the soliton itself. 

We have examined how the form of the localized solution (soliton) changes with the value 
of X2- This is illustrated in Fig. 1, where we see a gradual sharpening effect at the peak of the 
soliton solution starting from the standard Davydov form at \2 — 0. The soliton disappears 
at a relatively low value of \2 — —1.34pN, well within the range of physically acceptable 
values. This effect has been carefully investigated through analysis of the numerator and 
denominator of H(a) appearing in eq. ( 14.181) . As previously mentioned, a± = is required 
for solitary wave solutions to exist. Of utmost importance to what follows is the location 
and reality of the roots of the numerator in H(a). Obviously a double root exists at a = 0, 
and the remaining four roots are pairwise symmetrical with respect to sign reversal. We 
shall label the squares of these roots simply as a 2 =ri,r 2 . With this in place, we note that 
the singularity in the denominator of H(a) coincides with r\ = r 2 at a point \2 = X2 " ■ We 
have determined that the value of X2 & t which the solitonic solutions are destroyed is given 
by 

cr = 2Joxi (52) 

3A-^ + 2J 

Substituting the parameter values from eq. (5.1) into eq. (5.2) we find that for our system 
X c { = —1.34 pN which is close to what is physically expected. The role of X2 in changing 
the topology of H(a) is clearly illustrated in Fig. 2. We show in Fig. 3 the region of stable 
solitons as deduced from the relative position of the roots of the numerator of H(a) and the 
zero of the denominator. 



VI. CONCLUSION 

This paper has discussed the distance dependence of exciton hopping in the Davydov 
model. Two cases needed to be discussed separately. In the case of compressional solitons 
(Xi > 0), soliton stability was enhanced simply because in the soliton, the average distance 
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between molecules was decreased, leading to a higher hopping probability. If this was the 
case in peptide chains, the conclusions reached through thermal stability studies of the 
standard Davydov model (Ref. J3J and others therein) might be unduly pessimistic. 

The second case of xi < corresponds to a dilatational soliton. Here there exists 
a critical value of X21 beyond which Davydov-like solitons cease to exist. This critical 
value is surprisingly close to physical values and therefore it is to be expected that thermal 
fluctuations at physiological temperatures would destroy the soliton. 

By including this important physical effect, we have arrived at a new set of equations 
which possess solitons. These new solutions are clearly relevant to the physics of poly-peptide 
chains. 

ACKNOWLEDGMENTS 

This research was supported by an East- West collaboration grant from NATO, the 
Alexander von Humboldt Foundation and Deutscher Akademischer Austauschdienst. J.A.T. 
expresses his thanks for the hospitality and kindness of the staff and faculty members during 
his stay at the University of Warsaw. The authors wish to thank the referee for relevant 
criticism which led us to understand the role of the signs of Xi an d X2- 



11 



FIGURES 

FIG. 1. The plot of soliton profiles for various values of \2 (hi pN). The dashed line is the 
Davydov soliton (\2 = 0). Note how for \2 = —1.33 pN, the soliton becomes sharply peaked. 

FIG. 2. The change in the plot of H{a) for various values of \2 (hi pN) corresponding to those 
used in Figure 1. The plot with \2 = is the standard Davydov case. 

FIG. 3. The location of squares of the roots, T\ and ti, of the numerator in H(a) (see eq. 
(4.18)) and the singularity (dashed curve) in the denominator of H(a) plotted as a function of \2- 
Note that xf" delineates the boundary of the soliton region. 
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